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TIDAL ESTIMATION IN THE PACIFIC WITH APPLICATION TO SEASAT ALTIMETRY 


INTRODUCTION 

This is the second report of an investigation concerning the application of Proudman functions to the es- 
timation of tidal components in the world’s oceans. The first report presented results obtained in the Atlantic- 
Indian Oceans (Sanchez et al., 1986). 

This investigation focuses on the Pacific basin and incorporates some aspects which were' lacking in the 
previous investigation in the Atlantic-Indian Oceans. 

The application of satellite altimetry to the investigation of dynamical phenomena in the oceans is be- 
coming a well-defined engineering and scientific field. The ocean tides produced by the gravitational attrac- 
tion of the moon and the sun constitute an interesting and important aspect of these investigations. The dy- 
namical response of the oceans to the astronomical forcing functions is not an equilibrium one. This is easily 
understood in the context of the normal modes of the fluid basins involved, which have periodicities spanning 
the range of the forcing functions. This is certainly true for the gravitational normal modes modified by rota- 
tion. in relation to the diurnal and semidiurnal components of the tidal potential. The rotational modes modi- 
fied by gravity have longer periods; they might play a role in the response of the oceans to seasonal winds and 
other phenomena associated with ocean circulation. 

The frequency spectrum and space structure of the normal modes of the world’s oceans are determined 
by the configuration of their boundaries and their bathymetry, as well as by the physical properties of the 
fluid, such as density and viscosity. The actual detection of the free periods of oscillation of the oceans from 
the analysis of available data remains an open area for investigation. The reader is referred to Platzman, et al. 
(1981), and to Luther (1983), for excellent discussions of this problem, including further references. The 
numerical computation of normal modes for mathematical models of the ocean basins is an area where con- 
tributions have been made by many investigators. The following is not intended to be an exhaustive historical 
review, but only a sampling of the literature: the theoretical papers by Longuet-Higgins (1968), Longuet- 
Higgins and Pond (1970), and Christensen (1973), which dealt with the analytical eigensolutions of the Laplace 
Tidal Equations (LTE) for rotating spherical and hemispherical regions, Gotlib and Kagan (1980, 1982) 
reconstructed the spatial structure of the Kj and 0j tides from eight eigenfunctions of the LTE with periods 
in the range 17-30 h; they used 639 data points from coastal and island stations. Webb (1980) studied the role 
of resonance and friction in the tidal response and energy dissipation. He used a flat-bottomed, hemispherical 
ocean, centered at the Equator. Platzman (1978, 1984), and Platzman, et al. (1981), computed a range of 
normal modes for the world oceans and used them to synthesize some of the diurnal and semidiurnal tides. In 
a more limited geographical scale, Wubber and Krauss (1979) computed periods and structures of the gravita- 
tional free oscillations in the Baltic Sea; Schwab and Rao (1983) computed several of the lowest gravitational 
normal modes in the Mediterranean-Adriatic system. 

The objectives of this investigation are not limited to the numerical computation of normal modes; the 
main goal is the exploration of the applicability of the Proudman functions as a basis for the analysis of tidal 
data. However, the computation of the modes is still considered to be of sufficient scientific interest in itself 
to warrant the partial presentation of the results. In addition, the computation of the forced solutions and the 
potential energy spectra for the different tidal components contribute to the understanding of the method 
and lends unity to its presentation. The detailed mathematical formulation of the method can be found in 
previous publications (Sanchez et al., 1985, 1986). The theoretical foundation is Proudman’s theory (1918) as 
formulated by Rao (1966). The theory provides the formalism for calculation of the gravitational (first class) 
modes and the rotational (second class or Rossby wave) normal modes of irregularly-shaped basins with realis- 
tic bathymetry. 

The method requires the solution of two elliptic partial differential equations with second-order oper- 
ators, which are simpler than the tidal operator. The boundary conditions correspond to vanishing of the 
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stream function and normal derivative of the velocity potential. The elliptic operators are numerically rep- 
resented in finite difference form; the grid used is a Richardson lattice, which preserves self-adjointness. The 
solutions yield the velocity and surface height fields in terms of orthogonal functions with time-dependent 
coefficients. These functions are then substituted into Laplace’s tidal equations: if the homogeneous equa- 
tions are used, one obtains the normal modes; if the forcing terms are included, then the forced solution is 
obtained. In both cases, the solution is obtained numerically. The surface height field is only dependent on 
the velocity potential orthogonal functions, also referred to as Proudman functions. They correspond to the 
gravitational normal modes at zero rotation. The expansion coefficients of these functions can be estimated 
in a least-squares sense from available selected tidal measurements. It is possible to use the normal modes as 
a basis for the interpolation (Woodworth and Cartwright, 1986). We believe the adoption of the Proudman 
functions as a basis offers the advantage of simplicity, since to compute these functions, it is only neces- 
sary to solve the equation 


V * h V0 7 = - X 7 0 7 
subject to the boundary condition 



where 


h: 

V 

V 
is 

an 


variable depth of the fluid in equilibrium 
eigenfunctions of the velocity potential 
corresponding eigenvalue 

gradient in the direction normal to the boundary. 


Besides simplicity, the Proudman functions offer the advantage of being free of the dynamic prejudice 
introduced by factors such as friction. Also, they are dependent only on the shape and depth of the partic- 
lar basin, so they have to be computed only once when analyzing different tidal constituents. 


NORMAL MODES AND FORCED SOLUTIONS 

The numerical solutions were obtained by means of finite differences using a 6° x 6° Richardson lat- 
tice. There are 510 velocity potential points (4>) and 455 stream function points {$) in the Pacific basin. The 
eigenfunctions of the velocity potential, or Proudman functions, correspond to the gravitational normal 
modes of the basin at zero rotation. The periods of oscillation of the slowest modes are given in the first 
column of Table 1. These modes also have the largest space scales. The second column of Table 1 exhibits the 
periods of the slowest 30 gravitational modes modified by rotation. They all have kinetic: potential energy 
ratios close to unity, as shown in the last column. The normal mode solution was obtained by including the 
lowest 150 eigenfunctions from each field ( <t> and \li). Platzman (1978) computed normal modes for the world 
oceans using a finite-element technique with triangular grid elements having an average area equal to that of a 
4.5° equatorial square. A direct comparison with his results must consider the fact that our model has arti- 
ficial boundaries across the Drake Passage and the Tasman Plateau; a direct consequence of these barriers is 
the impossibility of obtaining any circulation around Antarctica. However, some similarities in the results of 
these separate investigations are clearly discernible. Platzman’s mode 16 with natural period 28.7 hours 
(Platzman, et al., 1981), which is a combination of a Pacific 1/2 wave and an Antarctic Kelvin wave, bears a 
strong spatial resemblance to our second gravitational mode with a period of 32.82 hours. The correspondence 
is verified still further by Platzman’s comments concerning the configuration of this mode at zero rotation and 
its period of 29.9 hours when using a model of the Pacific Ocean alone. Figures 1 a and 1 b show the Proudman 
function , and the corresponding normal mode obtained with our model. The Proudman function plot shows 
the contours of equal amplitude (normalized to maximum value = 100). The contour of zero amplitude is 
the nodal line; the nodal lines remain fixed in space and the nonrotating modes have the nature of standing 
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Table 1. Gravitational Modes 


Period (hr) 

Energy Ratio 

Nonrotational 

Rotational 

(Kinetic/Potential) 

43.91 

43.94 

0.97 

29.99 

32.82 

0.98 

26.37 

25.16 

1.01 

22.11 

21.61 

0.92 

20.63 

20.76 

1.01 

18.58 

18.65 

1.06 

17.93 

16.70 

0.99 

16.47 

16.12 

1.04 

16.10 

15.33 

1.00 

15.77 

15.00 

1.16 

14.53 

14.21 

1.33 

14.13 

13.77 

1.14 

13.18 

12.60 

1.00 

12.69 

12.26 

1.08 

12.46 

11.78 

1.23 

12.09 

11.57 

1.15 

11.66 

11.28 

1.17 

11.06 

10.84 

1.18 

10.92 

10.61 

1.07 

10.53 

10.24 

1.05 

10.09 

9.95 

1.05 

9.73 

9.38 

1.06 

9.51 

9.33 

1.17 

9.49 

9.10 

1.15 

9.19 

8.99 

1.12 

9.01 

8.79 

1.08 

8.80 

8.66 

1.10 

8.68 

8.39 

1.07 

8.43 

8.21 

1.08 

8.40 

8.09 

1.12 


waves. In the normal mode plot, the contours of equal amplitude and phase are given by the solid and dashed 
lines; the arrows indicate the sense of progression for high and low water. Platzman’s mode 19 with a period 
of 21.2 hours resembles our third gravitational mode with four cyclonic amphidromes and the period length- 
ened to 25.16 hours. Figures 2a and 2b refer to this mode. As the structures of the modes become finer, and 
their periods shorter, it becomes harder to establish clear cut one-to-one correspondences between the results 
of the two models. 

The forced solutions yield the amplitude and phase fields for the various tidal components. In the 
context of this investigation, the forced solutions are frictionless, without self-attraction or tidal loading 
effects. Figures 3 and 4 show the tidal solutions obtained for the M 2 and 0] components. It is interesting to 
note that the M 2 amphidromic pattern in the northeast Pacific resembles that obtained by Brown and Hutch- 
inson (1981) from analysis of SEASAT and GEOS-3 altimetry data, although the resemblance of both to 
Schwiderski’s ( 1 983) solution, Figure 5, is not as good. 

Tables 2 and 3 show the periods and potential energy percentages of the 10 most energetic modes for 
each of the diurnal and semidiurnal components. The mode with a period of 11.78 hours is the most en- 
ergetic for the M 2 and N 2 components. For the S 2 and K 2 solutions, the most energetic mode has a period 
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Table 2. Potential Energy Spectrum for the Semidiurnal Components 


• 

m 2 


s 2 


N 2 


K 2 


Rank 

Period(hr) 

P.E.(%) 

Period(hr) 

P.E.(%) 

Period(hr) 

P.E.(%) 

Period(hr) 

P.E.(%) 

1 

11.78 

19.76 

11.57 

27.47 

11.78 

43.5 

11.57 

24.18 

2 

12.26 

17.34 

11.28 

19.38 

13.77 

12.03 

11.28 

17.34 

3 

11.57 

14.33 

10.84 

13.04 

11.28 

9.15 

12.60 

15.45 

4 

11.28 

11.93 

12.60 

12.11 

9.38 

9.11 

10.84 

12.05 

5 

13.77 

5.47 

12.26 

7.32 

12.60 

4.64 

12.26 

9.14 

6 

10.61 

4.86 

11.78 

3.06 

16.12 

4.11 

11.78 

3.83 

7 

15.00 

4.86 

10.61 

3.05 

' 9.10 

2.07 

8.79 

2.55 

8 

10.24 

2.63 

8.79 

2.71 

10.84 

1.75 

10.61 

2.52 

9 

10.84 

1.80 

9.38 

2.09 

10.61 

1.71 

9.38 

1.78 

10 

16.12 

1.78 

8.99 

1.11 

8.99 

1.45 

13.77 

1.56 




Table 3. 

Potential Energy Spectrum for the Semidiurnal Components 



Kl 


0l 


P 1 

Qi 


Rank 

Period(hr) 

P.E.(%) 

Period(hr) 

P.E.(%) 

Period(hr) 

P.E.(%) Period( hr) 

P.E.(%) 

1 

16.70 

26.64 

25.16 

37.90 

16.70 

25.07 25.16 

31.82 

2 

21.61 

20.13 

20.76 

11.26 

21.61 

21.34 16.70 

15.04 

3 

20.76 

18.37 

16.70 

8.60 

20.76 

17.34 20.76 

12.07 

4 

16.12 

10.36 

18.65 

8.03 

16.12 

11.06 18.65 

6.52 

5 

14.21 

4.71 

15.00 

7.42 

14.21 

4.76 15.00 

5.96 

6 

13.77 

3.26 

12.60 

4.79 

43.94 

3.14 13.77 

5.50 

7 

12.26 

3.20 

13.77 

3.79 

12.26 

3.07 12.60 

5.22 

8 

43.94 

2.62 

21.61 

3.08 

13.77 

3.04 12.26 

4.81 

9 

18.65 

1.89 

12.26 

2.65 

25.16 

2.20 43.94 

3.39 

10 

25.16 

1.51 

43.94 

2.19 

18.65 

1.61 14.21 

2.89 


of 1 1.57 hours. The energy percentages given in the table are relative to the total (100 percent) contributed 
by the 1 50 modes in the solution. For the diurnal components, the most important modes are the ones with 
periods of 16.70 hours (Kj and Pj), and 25.16 hours (Oj and Qj). The rms amplitudes computed for the 
various components are as follows: M 2 (63.9 cm), N 2 (25.2 cm), S 2 (24.2 cm), K 2 (6.9 cm), Kj (10.9 cm), 
P] (3.6 cm), Oj (8.3 cm), Qj (1.1cm). 
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FITTING SCHWIDERSKI’S NUMERICAL SOLUTIONS 


Schwiderski’s numerical solutions to Laplace’s tidal equations, constrained to fit all the measured 
coastal data (Schwiderski et seq., 1983) provide a way of testing the efficiency of the Proudman function 
representation, vis-a-vis a spherical harmonic expansion approach. Schwiderski’s amplitudes and phases for 
the various components have been used as data at the 510 velocity potential points of the grid. Figures 5a 
and 5b show Schwiderski’s solutions for the M 2 and 0j tidal components. The tidal fields were then fit by 
least squares, first in terms of Proudman functions, then in terms of spherical harmonics. It is expected that 
the Proudman functions will yield better results when fitting diurnal and semidiurnal components which 
have periods in the range spanned by the gravitational modes. However, the expected outcome is not as clear 
when fitting tidal components with long periods. In that context, the lunar fortnightly (Mf). the lunar 
monthly (M m ), and the solar semiannual (S^) have been included in the tests. The results are given in Table 
4. The second column lists the rms for the amplitude of the various numerical solutions; the third and fourth 
columns give the results of the fits. The rms of fit in the last two columns of Table 4 is defined as follows: 


510 



(AEj - ASp 



Vi 


where AEj denotes the amplitude of the estimated solution at point i, and ASj stands for the amplitude from 
Schwiderski’s solution. As expected, the Proudman functions perform better (i.e., smaller rms) for the 
diurnal and semidiurnal components; the M 2 component can be estimated better with 150 Proudman func- 
tions than with a 14th degree spherical harmonic expansion involving 1.5 times the number of coefficients. 
The fortnightly and lunar monthly were estimated about equally by both methods, but the spherical har- 
monics give better results for the estimation of the solar semiannual component. In this last case, a 9th 
degree expansion (200 coefficients) yields better results than the fit using 450 Proudman functions. 

The potential energy spectra corresponding to the fits in Table 4 are partially given in Table 5 for the 
M 2 and 0] components, for the cases in which 150 Proudman functions were used. Comparing those with 
the results given in Tables 2 and 3 for the theoretical solutions, one sees some similarities and some differ- 
ences. Seven out of the first ten modes, and two out of the first five, are common to both M 2 fits. In the 
case of 0j, eight out of the first ten, and three out of the first five, are common to both fits. Figures 6a and 
6b show the solutions for M 2 and 0j when fitting Schwiderski’s solutions with 200 Proudman functions. 
Schwiderski’s numerical solutions were also fitted using a reduced data set, consisting of 362 points covering 
the central part of the Pacific Ocean (Figure 9a). These are the points used to analyze the SEASAT altimetry 
data as explained in the next sections. Schwiderski’s M 2 tidal map was fitted using the first 90 Proudman 
functions (in order of decreasing space scale and length of period); an rms of fit (over the 362 data points) 
equal to 3.10 cm was obtained. Schwiderski’s 0] map was fitted using the first 50 Proudman functions; the 
rms of fit obtained is equal to 0.6 cm. The coefficients derived from the fits above were used to compute the 
potential energy spectrum, and the results are given in Table 6. The potential energy percentages are based 
on a 100 percent total using 90 modes for M 2 and 50 modes for 0j. For comparison, the M 2 and 0j theo- 
retical solutions obtained by using 150 modes were also fitted. The corresponding rms of fit were 0.57 
cm for M 2 (using 90 functions) and 0.17 for 0j (using 50 functions); both over the 362-point data set. Tire 
corresponding potential energy spectra are also given in Table 6. The first thing to notice in these results 
is the preponderance of one mode in all cases— the mode with a rotational period of 21.61 hours. The cor- 
responding Proudman function has a nonrotational period of 22. 1 1 hours. This mode is especially impor- 
tant when fitting the 0] component; it then accounts for more than 50 percent of the potential energy. In 
general, the fits to Schwiderski’s solutions and to the theoretical solutions show similar patterns in the poten- 
tial energy spectrum although they are not identical. Six out of the first ten, and three out of the first five, 
are common in the fits to M 2 ; seven out of the first ten are common in the 0j fits. Also, the first three are 
common, and in the same order, in the fits to Oj. 
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Table 4. Fitting Schwiderski's Numerical Solutions Using Proudman Functions 
and Spherical Harmonics. RMS Values in cm. 


Tidal 

Component 

RMS Amplitude, 
Schwiderski 

Proudman functions 
RMS of Fit (no. 
of functions) 

Spherical Harmonics 
RMS of Fit (degree 
of expansion) 

m 2 

33.82 

4.40(150) 
3.54 (200) 
1.03 (450) 

8.81 (9th) 
5.56 (14th) 

s 2 

11.67 

1.29(200) 
0.66 (450) 

2.98 (9th) 
1.62 (14th) 

*1 

17.16 

2.48 (200) 
1.17 (450) 

4.40 (9th) 
2.91 (14th) 

Ol 

12.19 

2.01 (150) 
1.74(200) 
0.71 (450) 

3.17 (9th) 
2.14 (14th) 

M f 

1.21 

0.18(200) 
0.13 (450) 

0.22 (9th) 
0.15 (14th) 

M m 

0.90 

0.13 (200) 
0.05 (450) 

0.21 (9th) 
0.14 (14th) 

^sa 

1.52 

0.43 (200) 
0.39 (450) 

0.30 (9th) 
0.17 (14th) 


Table 5. Potential Energy Spectrum for Schwiderski's M 2 and 0 1 Using 

Fits to 510 Points 


Rank 

M? 

Period (hr) z 

P.E. (%) 

Oj 

Period (hr) 

P.E. (%) 

1 

11.57 

22.10 

18.65 

17.90 

2 

12.26 

12.87 

25.16 

15.66 

3 

12.60 

6.08 

16.70 

13.73 

4 

10.84 

4.45 

13.77 

8.29 

5 

10.61 

4.26 

32.82 

7.27 

6 

15.00 

4.08 

16.12 

5.73 

7 

16.70 

4.03 

12.26 

3.61 

8 

11.78 

2.53 

43.94 

2.74 

9 

15.33 

2.42 

12.60 

2.68 

10 

11.28 

2.20 

21.61 

2.47 
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as mapped in 







Table 6. Potential Energy Spectrum for M 2 and 0^ Using Fits to 362 Data Points 


M 2 (Schwiderski) M 2 (theoretical) 0j (Schwiderski) Oj (theoretical) 


Rank 

Period(hr) 

P.E.(%) 

Period(hr) 

P.E.(%) 

Period(hr) 

P.E.(%) 

Period(hr) 

P.E.(%) 

1 

21.61 

15.56 

21.61 

19.35 

21.61 

56.76 

21.61 

57.77 

2 

9.95 

11.36 

20.76 

14.60 

20.76 

20.75 

20.76 

18.80 

3 

11.78 

10.22 

9.95 

8.98 

25.16 

3.43 

25.16 

9.57 

4 

14.21 

8.24 

11.78 

6.54 

16.70 

2.76 

18.65 

3.69 

5 

9.33 

5.47 

18.65 

4.45 

12.60 

2.24 

16.12 

1.70 

6 

10.61 

4.58 

12.60 

3.95 

32.82 

1.87 

32.82 

0.97 

7 

25.16 

4.56 

10.84 

3.71 

18.65 

1.58 

15.00 

0.87 

8 

8.66 

3.59 

13.77 

3.39 

10.61 

1.52 

15.33 

0.64 

9 

20.76 

3.52 

10.61 

3.07 

9.10 

1.02 

11.28 

0.57 

10 

10.24 

2.89 

14.21 

2.70 

11.28 

0.97 

10.61 

0.55 


TIDAL ANALYSIS OF SEASAT ALTIMETRY 

Expansion in terms of a proper set of spatial basis functions is, in principle, the best way of extracting 
the tidal field from sea surface topography measured by a satellite altimeter. Attempts have been made to 
define the M 2 tide in the global ocean by analysis of the altimetry from SEASAT using expansions of spheri- 
cal harmonics (Mazzega, 1985) and of computed normal modes (Woodworth and Cartwright, 1986), both 
with some degree of success and with some limitations. The limitation of spherical harmonics is that, since 
data noise severely limits the order of expansion, it cannot properly model the finer tidal stnicture near the 
edges of the ocean basins, nor the abrupt transition from one ocean basin to another across a narrow land- 
mass (e.g., Central America). Computed normal modes of oscillation correctly model the requirements of 
ocean basin geometry, but are limited by the physical approximations (e.g., lack of friction) built into their 
formulation. Both Mazzega (1985) and Woodworth and Cartwright (1986) achieved better results by separat- 
ing the solutions for different ocean basins at the expense of discontinuity across their common boundaries. 

The set of Proudman functions for the Pacific Ocean described in the first sections of this paper should 
provide a better alternative to both spherical harmonics and normal modes, because they embody the correct 
requirements of mass conservation and boundary flow without imposing any artificial dynamic constraints 
in their formulation. An expansion of the tidal field in Proudman functions may, therefore, embody any 
reasonable regime of linear internal dissipation, which is important in tidal theory. The principal limitation 
is the implied lack of energy transmission across the model boundaries, but this is not considered to be 
very great in the seas bordering the Pacific Ocean. (The Bering Sea, at one time considered to be a major 
sink of M 2 tidal energy, has since been shown to be insignificant in this respect (Sundermann. 1977)] . 

In applying Proudman functions to analyze SEASAT altimetry, we restricted the analyisis to an in- 
terior zone of the Pacific (Figures 9a, b) which covers about three-quarters of the total area used for comput- 
ing the functions, in order not to confuse the fitting by fine tidal structure near the coastal boundaries. In 
particular, the difficult tidal area around New Zealand has been excluded, as well as regions of high meso- 
scale activity in the Kuroshio Current off Japan and at high southern latitudes, (Cheney, Marsh and Beck- 
ley, 1983). The loss of ideal orthogonality of the set of functions over the smaller area is not too serious, and 
in any case, the area of integration is weighted by the disposition of valid data along the satellite’s ground- 
tracks. The use of nonzero cross-products between members of the set is a normal part of the analysis pro- 
cedure. 

The Altimeter Data Set 

We used a file of altimetric Sea Surface Heights (SSH) and related data specially prepared from the 
original GDR tapes by Oscar L. Colombo of EG&G for purposes such as tidal analysis. All data were taken 
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from the ‘3-day repeat’ phase of the SEASAT mission (1978, Sept. 13 - Oct. 10), here reduced to the period 
Sept. 15 - Oct. 9, as far as roughly-continuous data series are concerned, apart from some gaps of a few 
hours. The pattern of ground-tracks for this period, repeated every 3 days, with a minimum separation of 
8.4° in longitude between parallel tracks, has been reproduced many times; see, for example, Figure 4 of 
Brown and Hutchinson (1981). 

The original Is - 1 data had been smoothed by Colombo to roughly 10s intervals corresponding to 
about 70 km of ground-track, and adjusted to a pre-set grid of earth-positions. The altimeter heights 
along tracks which, owing to a slight drift in the orbit, did not quite pass over the pre-set positions, were 
corrected by an adjustment in height proportional to the cross-track slope of the geoid in the region con- 
cerned. All data which, for a variety of technical reasons, showed gross departure from consecutive smooth- 
ness were eliminated from the data set. 

The computed orbit used to correct the altimeter heights to SSH was the Goddard Space Flight Cen- 
ter’s PGSS4 orbit, which was more accurate than any other computations for SEASAT prior to 1987. Apart 
from this, a standard suite of algorithms for the correction of radar path-length had been applied to the data 
at an early stage of processing. These corrections included subtraction of the ‘body-tide’ of the earth’s 
crust and a computed ‘ocean tide’ from the models of Schwiderski et seq. 1980. Since we are concerned 
with analysis of the ocean tide signal, we restored it by adding to the given SSH data exactly the same fig- 
ures for ‘ocean tide’ which had previously been subtracted. The geocentric tidal signal contained in the 
adjusted SSH consisted, therefore, of the ocean tide relative to earth as normally recorded, slightly modi- 
fied by the ‘load tide’, or ocean-induced deformation of the ocean floor, (Parke and Hendershott, 1980). 

Corrections for Geoid and Orbit Error 

The adjusted SSH, as described above, is dominated by the geoid signal of order ± 50 m and to a lesser, 
but important, extent by the error in the computed orbit, of order ± 1 m. These have to be subtracted by _ 
some means for oceanographic analysis. The repeated track provides the most convenient means of removing 
at least the geographically-invariant part of the SSH, consisting of the geoid, though somewhat distorted by 
the quasi-permament dynamic ocean topography. A commonly-employed method is to reduce each recorded 
height SSH (0, X, t) at time t to the anomaly 


DH ( 0 , X , t) = SSH ( 0 , X , t) - MSH ( 0 , X ) (4.1) 


where MSH (0, X) is the mean of all recorded heights along repeated tracks over the position with latitude 
0, longitude X. This is a good method in principle when applied to a very long record, but in the case of SEA- 
SAT, with < 8 repeated passes over each section of track, frequent irregular fallout of some data values intro- 
duces noise in MSH, owing to the irregular distribution of orbit errors among the eight possible passes. Such 
noise largely accounts for the irregularity of the histograms of anomalies shown in Figure 6 of Woodworth 
and Cartwright (1986). 

We have avoided such irregularity in the present study by working in terms of the ‘running difference’ 


AH ( 0 , X , t) = SSH ( 0 , X , t) - SSH ( 0 , X , t + KT), 


(4.2) 


where T closely approximates the mean repeat period, 3.0089 days in the case of SEASAT (the exact 
value being determined by exact repeatability of 0,X), and K is an arbitrary small integer. In fact, K was set 
equal to 2, because the time difference of 2T conveniently matches the 6d duration of the consecutive arc- 
lengths used for computation of the PGSS4 orbit. Whenever either of the two variables on the right of (4.2) 
was missing or invalidated, the corresponding value of AH was omitted from subsequent analysis. 

AH, as defined by (4.2), may not be a useful variable for more general oceanographic analysis, but it 
is well suited to tidal analysis. Omitting (0, X) for brevity, a harmonic tidal component 

SSH(t) = H cos wt 
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appears in AH as fH cos (wt - <t>), where 

f = 2 sin (w KT)/2, <t> = {it - « KT)/2 (4.3) 

Values of f and <p for six leading harmonic tidal constituents, with K = 1 , and 2 and T equal to the SEASAT 
repeat period, are listed in Table 7. The lunar constituents Qj, Oj, N 2 , and M 2 are enhanced in amplitude, 
while the solar constituents Kj and S 2 , being aliased into very low frequencies by the near-sun synchronism 
of the orbit, are conveniently suppressed. The corresponding lunar harmonic constituents in SSH(t) are, 
therefore, easily recovered from their values in AH(t) by applying well-defined adjustments to amplitude 
and phase. 

We still have to consider the removal of error in the computed satellite orbit, which in the case of AH 
with K = 2 is, in effect, the difference in orbit error between two consecutive 6d arcs. Orbit error is known 
to have very long spatial scales, and the usual approach to limited areas of ocean is to subtract an arbitrary 
‘bias and tilt’ from each pass, calculated to minimize DH(t) or some set of crossover differences. Such meth- 
ods are not appropriate to the analysis of long wavelength oceanic phenomena over large ocean basins, be- 
cause some of the signal will be eliminated as part of the bias and tilt. (Woodworth and Cartwright ( 1 980) 
describe methods of tidal analysis in which the small-scale tilt of DH is itself taken as the variable, thus par- 
tially removing the problem of bias.) We have tried a new approach in which a nearly sinusoidal form. 

E(t) = A(t)cos n t+ B(t) sinfl t + C(t) (4.4) 

is optimally fitted to the variable 

A 0 H(t) = AH(t) - AH'(t) (4.5) 


where A, B, C, are slowly ranging arbitrary functions of time, fi = 86 jt/T, the ‘orbital frequency’ for SEA- 
SAT. and AH' is the difference in (Schwiderski) model ocean tide corresponding to (4.2). 

It is well known that a large part of the orbit error is spectrally concentrated near the once-per-revo- 
lution frequency S2, but theoretical studies have shown that the spectrum is not an infinitesimally narrow 
line which would imply constant values of A, B, C for any arc or for each of a pair of arcs. It has also been 
suggested that A, B, and C could be represented with low-degree polynomials over ephemeris arcs up to 7 
days long; Colombo (1984) has proposed using quadrics for A and B and somewhat higher-degree curves 
for C. We have compromised by evaluating the mean values of A, B, C, over all valid global data AqH during 
sequences of three consecutive orbits; the origin of t taken centrally at an ascending equatorial crossing, and 
least-squares covariance matrices weighted by a cosine-taper window. A typical set of sequences of such 
values, evaluated with central origins at every consecutive orbit during the 6-day arc commencing 1978 
Sept. 1 7.0, is shown in Figure 7a. 

There is jitter in each sequence of coefficients, owing to the irregular sampling of ocean track along 
each orbit, to the spasmodic fallout of data, and to areas of residual tide from inaccurate AH'; but distinct 
secular trends are obvious, justifying the general form (4.4). For definition of the orbit error E(t) to be sub- 
tracted from every sample of AH(t), A, B, C were approximated by smooth polynomials up to the 4th 
order as depicted in Figure 7b. The curves are clearly not linear, and the amplitude (A- + B^) 1 / 2 depicted 
by a continuous line never reaches zero, but varies between about 0.5 and 1.0 meters. These figures accord 
with simple observations of cross-over differences in SSH from the PGSS4 orbit. 

Typical variances of the global data for AQH(t), before and after subtraction of E(t) from (4.4), are 
0.142 and 0.05 lm^ respectively, both figures being averages over the same several samples. Note that the 
second, i.e., residual, variance is expected to be about twice the intrinsic variance of the sea surface (tide- 
reduced) because of the definition of AH (4.2). Within our interior Pacific area, the corresponding figures 
from 19,351 values of AgH(t) were 0.039 and 0.008m2, respectively. These figures compare satisfactorily 
with the typical rms residual sea surface, about 0.07m, deduced for the same area from the map of Cheney 
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Figure 7a. Estimates of the coefficients A, B, C, in meters (equations 4.2, 4.4, 4.5), defining 1 rev”^ harmonic components of 
tide-reduced differenced sea surface topography derived from SEASAT, reflecting the principal component of the 
orbit error. Individual points are evaluated at the ascending node of each revolution by summation over 3 consecu- 
tive revolutions. Revolution numbers (abscissa) are arbitrary, with 31 corresponding to 1978 Sep. 17, 2.079h. 



Figure 7b. Low-order polynomial splines fit to the harmonic 1 rev 1 coefficients plotted in Figure 7a. The line labelled rms 
is the corresponding harmonic amplitude (A^+B^) 1 ^. Notation and scale as in Figure 7a. 
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and Marsh (1983). The variance figures for this area after restoration of the tidal signal, that is, for AH and 
for AH - E, were 0.235 and 0.202m2. The tidal signal is clearly the dominant variation in our corrected data 
set. 


Choice of Tidal Time-Harmonics 

Sampling at a time interval T imposed by the satellite’s orbit aliases the harmonic constituents of the 
tide, as sampled at any particular point on the ground track, into much longer periods T a , (see, for exam- 
ple, Parke and Rao, 1983). The values of T a in days for the dominant constituents are listed in the last 
column of Table 7. The duration of the repeated-track data, further reduced by effectively 6 days by the 
use of the variable AH (4.2), limits our data span available from analysis to about 17.5d. From Table 7, this 
covers a little more than one complete cycle of M2 or Oj and a little less than two complete cycles of N2 
and Qj. Only a small fraction of a cycle of S2 and Kj is covered, but these and other solar constituents may 
be ignored due to their small f factors. 

One cycle of oscillation makes a rather poor sample for harmonic analysis at an individual point, but 
the use of spatial basis functions greatly improves the sampling efficiency by relating data along and between 
tracks. The spatial sequence of tracks presents additional cycles of phase change. By a similar process, the 
diurnal and semidiurnal tides, which have rather similar values of T a , are decoupled by including both ascend- 
ing and descending tracks in the same spatial scheme. (The correction E(t) (4.4) is also very relevant to the 
latter case). 

In order not to overburden our analysis scheme with too many unknown parameters, we linked N2 and 
M2 as an anomalistic monthly variation of arbitrary amplitude and phase, and linked Qi and 0 j in a similar 
fashion. In exact terms, instead of seeking a solution for M2 in the form 

X) coscj t + B n (6 , X) sinw t) 

we replace cosa> t, sinw t by 


c(t) = a(t) coscj t + 0(t) sinw t, s(t) = a(t) sina> t - 0(t) cosu> t (4.6) 

respectively, where 

a(t) = 1 + e cos (6 t — 7) 

0(t) = e sin (5 t - 7) 

5 = 2jt/ 27.555 d” 1 (4.7) 

e, 7 are assigned constants. 

The elliptical modulation corresponding to the ‘equilibrium tide’ would give e = 0.191, 7 = 0 for both N2 
relative to M2 and Qi relative to 0i. A simple analysis of the relationships between these constituents in 
Schwiderski’s model for the area of the Pacific Ocean gave 

e - 0.20, 7 = 22° for M 2 , 
e = 0.21, 7=8° for Oj, 

both results having coherence greater than 0.9. We used these numerical values with (4.6, 4.7) in all the fol- 
lowing analyses. 

We have thus reduced the number of tidal constituents which are relevant to this data set to one semi- 
diurnal and one diurnal harmonic with independent parameters to be assigned by least-squares analysis. All 
previous tidal analyses of SEASAT altimetry have ignored the diurnal tides. We have retained them because 
the diurnal tides in the Pacific are larger than elsewhere with respect to the semidiurnal tides, and also be- 
cause our procedures for removing the geoid and orbit errors are much better than in any other published 
exercise of this kind. 
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Table 7. Magnifications and Phase Lags (Equations 4.2, 4.3) 


Darwin 

w/2jt 


K = 1 

K = 

2 

T 

*a 

Symbol 

cycles d~ * 

f 

<t>° 

f 

0° 

Days 

Qi 

0.8932 

1.6624 

-33.78 

-1.8484 

202.45 

9.6 

0l 

0.9295 

1.1916 

-53.43 

-1.9140 

163.14 

14.8 

K i 

1.0027 

-0.1074 

-93.08 

0.2144 

83.85 

176.4 

n 2 

1.8960 

-1.6004 

143.15 

-1.9197 

196.30 

10.1 

m 2 

1.9323 

-1.1036 

123.49 

-1.8408 

156.98 

16.4 

s 2 

2.0000 

0.1112 

86.81 

0.2221 

83.62 

171.1 


Expansion in Proudman Functions 

Formally, we evaluate the set of coefficients aj n , bj n , b^n. n = 0(1) N — 1, which minimizes 
the mean square value of 

N— 1 
n = 0 

over the space and time parameters (0 , X, t) determined by the altimetry. Here, cj and sj are the modulated 
harmonic functions defined by (4.6), with constants particular to tidal species i, and P n , the computed 
Proudman functions of order n, including the ‘constant’ function of order zero and other values n increas- 
ing with ‘non-rotating frequency’ of the eigenfunction (Sanchez et al„ 1986). 

By standard theory (4.8) gives rise to a set of normal linear equations for the 4N parameters aj n , bj n 
which are solved by inverting a non-singular matrix of order 4N x 4N. A sequence of trial values N were tried 
up to N = 100, giving rise to pairs of variances v-(N), each being the contributions of the terms pertaining 
to tidal species i to the total variance V of AH, when N orders of Proudman functions are used. The residual 
variance of (4.8) is V — vj — Vj, and this, in general, decreases monotonically with N toward an asymptote 
representing the residual noise in the data AH. 

Results 



Skn C i< t > + b i,n s i^> 


P n (* 


X) 


i=H 


— AH ( 0 ,X ,t) 


(4.8) 


Figure 8a shows the variance v-,, split into partial components v c 2 , v s 2 associated with the functions 
c->(t), S 2 (t) respectively, from solution of the complete set of normal equations from (4.8), with N increasing 
from 5 to 100 in steps of 5. The ratio v^.'v^ depends on the time origin and is of no interest per se, but 

v c2 + v s2 = v 2 

would be independent of time and is the variance (scaled in cm2) contributed by semidiurnal terms. The 
total data variance V = 2020 cm2. 


As is usual in such least-squares fitting, the partial variances v increase rapidly for small values of N, 
then increase more slowly. Occasional slight decreases of v with N are due to numerical deficiencies in the 
matrix inversion procedure. It is hard to decide exactly for which value of N, further increase is unprofit- 
able; ultimately, the addition of higher-order functions merely tries to fit noise in the data with increas- 
ing ‘wiggliness’. A sign of over-fitting in this way is the appearance of large numerical values in the list of 
coefficients aj n , bj n ., Taking this factor into account, as well as the values in Figure 8a, resulted in the 
choice Of N = $0 as the largest reasonable compromise. 
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VARIANCES 







The corresponding values v c j, v $ j, for the diurnal tide are depicted in Figure 8b, with scale increased 
by 10. These vj components are appropriately much smaller than the V 2 components, and numerical noise 
is more evident at the larger scale. The large irregularities are clearly associated with similar irregularities 
in the v-» variances, and are evidently compensating for them. However, the same sort of general increase of 
vj with N is apparent, and from similar reasoning, the choice N = 50 is a sensible compromise, giving 

vj - 99, V 2 = 1576, V — Vj — V 2 = 345 cm 2 . 

If we had chosen N = 90, which was also considered, the corresponding variances would be 106, 1598, and 
316 cm 2 , respectively. 

It is interesting that the diurnal variance vj is considerably less than the residual noise variance 
V — vj — V 2 * This is due to the space-time filtering implicit in the normal equations, and is analogous to the 
extraction of a small but significant peak from a spectrum of noisy data. The energy of a spectral peak may 
be much less than the total energy across the whole spectrum. In other words, the residual noise in this case 
is concentrated in much higher wave numbers and frequencies than can be accommodated by Proudman 
functions of order up to n = 100 or thereabouts, with tidal periodicities. 

The standard tidal amplitudes and (Greenwich) phase lags Hj * (0, X), Gj (0, X) are given by 

Hi = (A 2 i + B 2 i ) 1/2 /f i , ^ = arctan (Bj/Aj) - *; + <t>\ (4.9) 


N— 1 

where (Aj, Bi)= ^ 
n=0 

and fj, <t > j are the adjustments for conversion from AH to SSH (4.2, 4.3) with a small additional correction 
for the 18.6y modulation of the tide potential. <2>'j is the phase of the harmonic potential concerned, 0j or 
M^, at the time t = 0 used in the computations. 

Contour plots of H:, Gj computed by (4.9) from the solutions for N = 50 are depicted in Figures 9a 
and 9b. Comparison of Figure 9a with the Schwiderski model M 2 , (Figure 5a), shows a high degree of sim- 
ilarity; the principal features being two large zones of maximum amplitude near the Equator and two zones 
of even greater amplitude coastal maximum in the Gulfs of Panama and Alaska, interspersed with four cen- 
tral amphidromes and a nearly amphidromic minimum in the north near 190°E. Figure 9a also reproduces 
the amphidrome east of Papua, New Guinea, but transforms the region of low amplitude southwest of 
Japan in Figure 5a into a true amphidrome. Schwiderski’s map suggests that the two central maxima have 
about equal amplitude of 55 cm. The SEASAT/Proudman function analysis gives the maximum near 230°E 
larger amplitude. The differences between Figures 5a and 9a are quite minor and neither can claim perfect 
accuracy anyway. 

The satellite solution for 0j (Figure 9b) again agrees in broad features with the Schwiderski solution 
(Figure 5b). This is even more remarkable, considering the generally lower amplitudes everywhere. We will 
not enlarge on the few points of dissimilarity, which are trivial. The fact that we can obtain any reasonable 
approximation to the 0j map with the observed noise background is an excellent indication of the methods 
of analysis used, and promises well for future results from the much longer series of altimeter data now be- 
coming available. 

Comparison of the distribution of potential energy among the Proudman mode numbers n deduced 
from the above solutions with those presented in Table 6 is less satisfactory. Mode n = 4, with rotational 
period 21.61h, again has the largest contribution to both M 2 and 0j, and some of the other modes listed 
in Table 6 also appear in the top ten here, but with different ranks. Other modes appear in the present solu- 
tion which do not appear in Table 6, and the distribution varies rather erratically with different values of 
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Figure 9a. M- tidal synthesis corresponding to 50 Proudman functions derived from altimeter height differences (4.2) over 
the reduced area of the Pacific Ocean shown. Solid lines show conventional tidal amplitude in cm; broken lines 
show Greenwich phase lags at 90° intervals. Arrows are attached to 0° phase contours and indicate the sense of 
wave propagation around the amphidromes. H: local maxima; L: local minima (amphidromic if phase contours 
intersect). 
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N in (4.8). The nonorthogonality of the Proudman functions over the reduced Pacific area used here, to- 
gether with the noise background, must cause aliasing between different modes, making the distribution 
of potential energy with n being rather meaningless in this case. 


Expansion in Spherical Harmonics 


In order to assess the relative merits of expanding the tidal fields in Proudman functions and in spheri- 
cal harmonics, we also analyzed the same altimetry data described in (4.1) through (4.3) by solving the 
normal equations for a' n , b' n , which minimize the mean square value of 



CiW + b/ 



Q n (fl.,X ) - AH(0 ,X , t) 


(4.10) 


where the modulated harmonic functions Cj, sj are the same as used in (4.8). Q n (6 , X) represents an ordered 
set of the real and imaginary parts of Associated Legendre Polynomials of degree 1 , order m, 


Qn 


Y £m(0 ’X ) 


21 + 1 (g-m) ! 

4rr (C + m )! 


Pg m (sinfl) e imX 


(4.11) 


PCm( x ) being the common Legendre polynomial, and the sequence of (1, m) corresponding to n = 0, 1 , 2,... 
being 


( 0 , 0 ) ( 1 , 0 ), R( 1,1) I(U), ( 2 , 0 ) R(2,l) 1(2,1),... 

with R, I denoting real or imaginary parts. Complete coverage up to degree and order k - 1 requires N = k-. 
Thus, N = 100 covers harmonics up to degree and order 9, while N = 50 covers degree and order 6 with the 
addition of Y 7 Q. 

The coefficients a', b' in (4.10) were solved in a very similar manner to a, b in (4.8) using the same se- 
quence of N. Results for the partial variances v'j, v'^ are qualitatively similar to those for the analogous 
vj, V 2 , except that numerical variation with increasing N is distinctly greater; probably due to greater correla- 
tion between spherical harmonics over the area considered— very much less than the whole sphere. However, 
v'j and v '2 are, on the average, slightly less, as in the typical comparisons listed in Table 8 . 

Table 8. Comparison of Predicted Variances (cm^) 


Proudman Functions Spherical Harmonics 


N 

V 1 

v 2 

v l+v 2 

v 'l 

v '2 

v', + v '2 

20 

94 

1449 

1543 

88 

1275 

1343 

40 

93 

1570 

1663 

93 

1537 

1630 

60 

101 

1582 

1683 

101 

1560 

1661 

80 

103 

1579 

1682 

101 

1570 

1671 

100 

100 

1579 

1679 

98 

1584 

1682 


Whether the slight improvements in predicted variance achieved by fitting the same number N of Proud- 
man functions are significant would be hard to say. Significant improvement must be felt mainly in the im- 
proved matching of the Proudman functions to the coastal configuration and bathymetry, but our reduced 
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area excludes the regions where these factors are important for other reasons. Maps of H, G computed from 
the coefficients a', b' look too similar to Figures 9a and 9b to be worth reproducing here, but there are quali- 
tative signs that the advantage afforded by the Proudman function expansion is indeed in the sea areas close 
to the boundaries of the chosen area. Further experiments with longer altimeter data sets in other ocean 
areas will be needed before one can finally decide whether the advantages of Proudman functions outweigh 
the computational effort required to produce them. 

Conclusions 

This paper brings together two areas of research; one consisting of the computation of Proudman func- 
tions and their associated normal modes in large ocean basins (Sanchez et al., 1986), the other exploring 
methods of extraction of ocean tides from satellite altimetric data (Woodworth and Cartwright, 1986). 

The Pacific Ocean, considered here as a closed basin, has the usual wide spectrum of eigenfrequen- 
cies, with gravitational natural periods from 43. 9h to 3.5h. Some of the normal modes spatially resemble 
those computed by Platzman, et al. (1981) for a global ocean, although closure of the connected region en- 
circling Antarctica prevents too close a comparison in frequency. 

As in other approaches to tidal modeling, the diumal and semidiurnal maps resulting from direct dynam- 
ic solution are only qualitatively good, partly due to the lack of friction and other unmodeled physical fac- 
tors. Much closer agreement with the known features of the ocean tidal maps is obtained when the Proud- 
man functions are constrained to a least-squares fit to empirical data, as represented by Schwiderski’s models. 
It might be suggested that the degree of fit with 200 Proudman functions shown in Figures 6a and 6b could 
be matched with any set of orthogonal spatial functions, but Table 4 shows that the approximation is much 
better than with an equivalent number of spherical harmonics. This result clearly reflects the advantages of 
using the kinematic constraints embodied in the Proudman functions, which are totally lacking in spherical 
harmonics. 

Modeling the tidal signal in the altimetry by Proudman functions is more limited by the noise level of the 
data, but the maps obtained with only 50 functions, (Figures 9a and 9b), are encouragingly close to the 
Schwiderski maps based on tide-gauge data, and are marginally but significantly better than the fit with the 
same number of spherical harmonics. It is particularly encouraging to see a plausible map for the low ampli- 
tude Oj component emerge from barely 17 days of SEASAT altimetry. Several months of altimetric data 
were previously thought to be necessary to separate the aliased diumal and semidiurnal constituents. The 
present result must be at least partially due to the improved methods of geoidal and orbital noise reduction 
employed here. In general, these results promise well for future applications of similar methods to altimetric 
data of much longer duration (GEOSAT, ERS-1, TOPEX-POSEIDON), using Proudman functions computed 
for the global ocean. 
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is, in all cases, better than the corresponding value for the same number of spherical harmonics. 

Before fitting the Proudman functions to the altimetry from the 3 -day repeat cycle of SEASAT, the data are 
processed by novel methods. The geoid component is eliminated by taking collinear differences at a fixed time-lag 
of 2 repeat cycles. Orbit errors are reduced by extracting the 1 rev“* component at every ascending node; this com- 
ponent varies slowly and nonlinearly in time. The spatial fitting process includes M 2 and 0| frequencies, both of 
which emerge with significant and realistic tidal mapping, but residual noise in the data limits the number of Proud- 
man functions to about 50-60 before showing signs of ‘over-fitting”. Fitting the same data by spherical harmonics 
gives marginally lower predicted variance for the same number of parameters. 
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